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ABSTRACT 


A Fortran IV program is developed for computing the percentage points of the Chi-square distribu: 
tion, The Newton-Raphson procedure is used with good initial approximations. The maximum number of 
iterations is four and in most cases only one or two are needed. An extensive table is given of the percentage 
points x2 in terms of va positive number, and Ql? ), the probability integral of the Chi-square distribu- 
tion. 


TABLE OF CONTENTS 


Foreword. 0... ee Pea iead wo eR hee wrates Biheg seep ah eaahcat nad wang ee aceon d rie yc a acts : 
AUSINICE 2 eS ew hehe tated se ae ua ins Cant S436 badd (6G, nae Boo GI, Hohe aoe ante oe hoe 
b Enrodictivitis oo 36 se. FSeeeen es Pt araaltnvety eas Ge abd and eke ere goth 
2. Initial Approximations 6.0 ..00.0.0. re ae Hams eee 4 ee ee eee Meer 
3. Criteria for Terminating (V- AR, pe SE Suite oe Sle eS PN ae 
4. Computer Program... .. Sea alls SOR ais Tae Seer E Ar ita base Sheet fangien stan : 
Blow Chart... 20... Rabi oleh eh ust ar ea bles cs eaeos eta st wae gee sentra ha ters er 
Se Ears in, Warters Table cor? og Sek ay eee cae tech Ba oe ea 2d GRA REE 
RUlCienees. <n tage ne eee eee Be Bee aes sth devia dy ee Bate wee eae Sed eeahty Go Batedya. Back Ais 
Appendices 
A. Table A Percentage Points of the Chi-square Distibution 00.000 ee 
B. Table Bo Percentage Points of the Chrsquare Distubution . 0.000000 000008. ese 
C. Program Listing... 00.0.0 ....0.0.000- ; ave S Os ae RS ee Gada ante : y 
Distoliiiont Bstis cas ices keto wee Dea hat be toys peng Goes ' eh witout at 
it 
eae er ee” a. ee heel ae SE Cay & eee Fs, . a! 


1. INFRODUCTION 


The primary objective of this paper is to descr be the effective use of the Newton-Raphson procedure 
(N-R) for finding the percentage points of the Chi-square distribution, For a given value of Q(vly-) and one 
. . >. is 
of v > 0, the value of yo is found, where 


Oly? y= (2"2 bw! fet 2-b ar, p> 0 ) 
+ 
2 


0<GObKAS 1. Aehy=b-OWh2)  QUl0)= 1. 2) 
Here F(+) denotes the complete gamma function. 


The tinal results of this study area Fortran FV program which ts fisted in Appendix C, and an eaten- 
sive table of percentage points with six correct significant digits. The program was used on the CDC 6700 
to construct the table which contains 72,000 entries. The entire calculation took only 7.5 minutes with an 
average of (N+ A) iterations per entry. The table is given in Appendix A. A very condensed version ts also 
given in Appendiy B primarily to exhibit a corrected copy of the table m1, pp. 984-985], The program 
was also employed to check a table of 5290 entries by Harter, [4. pp, 137-foL}. No errors larger than one 
unitin the last place were found. They are listedin Section 3. 


Because the notation is more convenient, we shall work with the incomplete gamma function nities, 
: . ‘ Ba ? 
Ala, x) and Gta.x). [2] instead of the corresponding Meh? band Gtely*), where 


ie a 
Musb® | eters Meee 43) 
Ita) 
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| — =f =! 2 Fh ». 
(ha, 8) ca ote tes EMD, ry) 
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The relitronsheps between @E) @2) and (3), (4) are well known and easy to show wang (9) 


Bor easy reference. we note that the GY - Abb acters te an iterative pracedtire for computing a rat y 
ofa tanction fib) Co The dt & Lise retate as peven by 


beet ty ALTE). Heh. iW) 


met 
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A Depecapee 


If convergence occurs so that 


lim ¢, =) fiv) = 0, (7) 


nwo 
and moreover if /’(") #0, then 


tim lings ale = QO O, ved (6, pp. 329.330]. (3) 


| 

2 
N00 

A glance at (8) shows that for /,,,, sufficiently close toy the erior in f,,4, 18 approximately proportional to 

the square of the error ings, Taig ts the attractive feature of the (V = A), and it applies in our case since 

f= dPidx. or dQ/dx. ty always positive, or negative, respectively for x > OQ, 1.¢.. never zero. 


In order to successfully use the GV - A) to find x trom Pa, x) or Ga, x), the following questions 
required affirmative answers: 


1) Can Aa, x) and Gta, x) be computed efficiently to high relative accuracy? 
2) Can a “good” initial approximation to x, denoted hereafter by xy. always be found? 
3) Will the CV - R) always converge tor such xy? 


Treating the third question first, itis eas to show from (3) that for fixed a & 1. Aa, x) ts concave tor 
ally > OQ. and that fora > 1 /¥@, x) has an intlection pomt atv =a > Lwith P convex tor OS x sta ~ Tand 


concave forx >a - 1. Similar properties hold for Q(a, x). Typical curves are shown in Figure 1, 


It can be shown acalvtically or tram a geometric interpretation of the (NV - RK) that convergence wall 
occud when & |b provided vp Sv and also wheng > Tix = @ =. Thus global convergence can always 
be achieved without mach difficulty. 


Nevertheless, experimen tition shows that if vy ts chosen simply to assure comvergence in the Ways it 
descrthed, the rat of coavergence in most cases will he stow. Therefore for better ofticnenes, ettorts were 
nude te detecsnne aa vg “close” tay even though the sufficient conditions noted above might not be 


Figure 1. Typical curves of Ata, x) and (Hu, x) fur fixed a. 
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satisfied, Convergence then would be established from computer results. This approach lead to a set of 
inital approximations depending on 2, Q, anda. such that fore 7 Pe 1- € (e = 5X 1073) no more than 
four iterations of the (NV - R) were needed to obtain v correctly to about 8 significant digits in every case 
tried including the 72,000 entries of the table in Appendix A. 


Kor the GV - 2) to work so well not only are “ood” starting values vy required, but also 2, Q and 
then first derivatives must he available te high precision. This presented no problem because of a previous 
study of ours, [2]. from which a computer program was developed giving 2 and @ to 12 signiticant digits 
toralla 0.4 °° 0. The high precision requirements are necessary, of course, When vos desired to many 
digits. In addition, they are also needed when # ts very close to zero or one otherwise the dast term on the 
right hand side of (6) could carry no reliable digits with disastrous ettects on the iterations. 


Clearly then. the three questions posed above can be answered in the affirmative, and consequently a 
computer program based on the (NV - 2) to find the percentage points of the chi-square distribution can be 
constructed. 


In the neat section we obtain the intial approximations for vw. in Section 3, the enteria for terminat- 
ing the (N= Adare specified. In Section 4, the computer program is described, and a flow ehart ty given, 


2. INERTIAL APPROXIMATIONS 


Out objective in thy scetion iy to give approximations tory which result in rapid convergence at the 
iN - Aytoa tinal estimate for, sav it. By rapid convergence we mean that no more than 4 iterations are 
necessary for about S correct signiticant digits i a*, 


We remark here that some of the analysis below is heuristie with the emphasis on plausibility, and 
attane etfttcenes inthe use at the (V+ AR) at expense of rigor, The usefulness of the initial approxima: 
tons is evaluated primarily by computer expermmentation, No effort p made to obriin assacuited error 
bounds on iy, = vl where vy denates the initial appronmiatiop, 


The input values ot Ata. x) and Gta wdare denoted by 2 and Q. respectively. Por the remainder ot 
this section, as wells Sections Jand 4, the flow chart an page EL will be helpful. 


Pee sete ab ital approumations fora are used depending on whether a bora baa 1, then 
yo =n QQ, 
Considering ab tiest, fet 
Be: OVW). (9) 


Tien tor doe 38, we fave eed bar aig ay assumed that ya > bo then witha < | 
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This inequality implies .v < .7985 contrary to the assumption x > 1. Thus with x <1 an approximation fo: 
Xp. say u, can be obtained from 


at aN ax? 
i laerperare’ (lear ister uereerenrcr tard (mS Oe] UEC 10 
Na. x) Mat ;| (a+ pit (a + 2)2! | pp. 6 | (10) 


by dropping all terms but the first, 1.2,, 
w= [Pla + Dy), (it) 
This result ts improved by using the first two terms of (OQ) with (11). ba this case 
tla 


my au oe te 
ocr Xo l- Pel = XQ l- ae 
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[sy w/t ula + i], B24. aS. | mies 


We sve that | - aa + PbO. Indeed, tram @10), <x forx <b. Hence ala + bp<ox/at b <x and 
b-uitat Lp > ben > 0. 


WAS 45 (witha <1), then an approximation tor vy is toand from the aay miptotic series 
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where 
Riad eT Se /V ad. (td) 
We biave by carving only the fist tenn on the mgt hand side at (13) 
eh yh Onan Be AS, | tS) 
anid taking bogarithins 
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Vora 21ND the first term on the lett dominates the second, giving as a ast approximation to vy 
ve ela kh. city 


‘Vins estimate tery, is then nproved by sna tin (16) to get 


vepe(a~ lee. (is) 


We improve the result sill further by using a rational ayoroxmmation for @, (56, 
, E Mi ( 
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and obtain 
v™ -dn Beta. ina ode flee 2 - av t tif (20) 
Now using boas given by € £5), for on the right hand side of (20), we fave 


Ip “y -th-aWn-in lis +2 -aniv# Ooi]. Boas, USES. (2) 


We digress here brietly to show that the arguments of the logarithms in (2) are always positive tor 
tumissibly values ofa. iv. Indeed, stnee B45, 
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40) TOSS tor alls 45.a< bo Abo. @ 2-0 7085 @ 2a > | 7985, 


itis possthle when wis it the neihborkoudt of the computer deterstned valine (Lol that ett 
better than él} dase hone tore, even thepuh Ae SS) or thre bbis better than 2b when Aes 48 
Bereding vy trom gia bs wy aid vy, from é 2b) by Nya twas conchuded trom commater expe rane titation 
that when Bo AS aod uy + ST he abs and Np? © 67 that usally the Best Ghote tor thes: sittataats 
wars te take vy ay the Larger Of egy ainda. 
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where the coctficrents are determined by the Chebyshes minimis error citenan  Phev ate give és 
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Computer results indieate that the mame renitive errer pio dareer than POS par dae Pe ae SY 
sitter expression is eiver tn Pho p O84} wath that athena fame tron con tinieg ome tess terme both 
numerator and denominator Per salnes of ob SR Ts wes feund tor bl - Ooand thers roghy ed hs 
fyb One vy computed, the Cornisirbisher espanston [3p vickds an approsination uy toe tarieh 
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The root in (31) with the plus sign, x(+) should not be used, because it does not yield a good 
approximation to x when x/(@ + 1) is small. Indeed, the approximation for w given in (29) is only good if 
x/a +1 is small. Therefore, we have from (29) ¥ =u. In addition, we also conclude from (30) that in fact 
u< X,so that 


ul(at+1)<Xhat t)<<1. (33) 


~ A 


Then expanding the radical of (31) gives the desired result, x 
precludes the use of x(+). 


u from x(-) but not from x(+) which 


The choice between Xg and Xq (x(~)) was made by experimentation. The conclusions were: If 
as 1/2.aS 19, Xo <O.11 (@+ 1) and a> 4, or instead of a > 4,a<4 and P< nj, where j equals that (= 1, 


ts 4 2.3.4 for which ia ~ ff is a minimum and ny = 15, my = .05,n3 = O1,n4 = .001, then Xq from (32) is taken 
sy as Xq. On the other hand xq of (26) is taken as xy under the conditions specified earlier, or when P< 1/2, 
5 Si aS 19,9 SO11(a + 1),a<4 and P>n,. 
aan 
Le = When P> 1/2 anda > 11, then again Xq isa satisfactory choice for xg. If however P> 1/2 and 
are . 1<a@< 11, then experimentation shows a better estimate than x9 is obtained from 
Aa 


: : E Xq =V-Cl- an NXg - ln [Gio > 2 -a)/(Xq + 1) (34) 


provided that in addition we also have 


tt = 
ze ne ; 2 
Le ay “ 
a Equation (34) is obtained in much the same way as (21) so its derivation is not given. 


3. CRITERIA FOR TERMINATING (N-2) 


in thas section We discuss the wavs in which the (V - A) is terminated. Since we are interested in 
_ SeUONEY corfect fo about & significant digits we required 


Ny Ryall axy, | ABSK 1G, (36) 


whete wg denotes the VM perate of (VR) WH the inequality holds, then xy ad xq y agree to atleast 4 
nig eaticant degts. However. this mean Ay has close ta Be since by (8) 


saith 
ivy > XlW DD Abeay + 80 tay? AS) we dex (37) 


« 


owhteh gngdics ihe eedative erron in yy is proportional Go the square al the relative error in Nyy. The 

" opertranality factor A can he levee when g and 7? ase large, for example, A 3525 wing : Woe ge torte, 
but sabes tee fer bige gd that yy every chew toy. En fact far the example above ay. tran (20), already 
kes Vito TE seuticant digits. Counputer results support aur conchision drawn feo ¢37). 


Although the primary goal is to insure the accuracy of x* to about ¥ significant digits (it is recalled 
X* 1s the final iterate from (V - R)). See page 3. A secondary objective is to satisfy (40; or (41) below. 


; * In this connection, we required once (36) was satisfied, that 
2 atl 4 ae 
wy - Nyy vy - Ma tl -ay) <2 _ ANS: ifa>xy-] (38) 
OT 
IX: - Nyy 7(xy + | . aX: + 2 = a) < 2AX A (YY + 1 ), if ad < Vay ~ 1. (39) 


If (36) 1s satisfied but the appropriate inequality, (38) or (39) is not. then one more iteration of (N= R) 
is carried out before stopping the procedure. We note in this case x* = x4, ALA= 5X 1075 no cases 
were found where the designated inequality, (38) or (39), was not satisfied, however for smaller values of 
A. for example 6 X 1077 some cases did occur and one more iteration was made. 


Inequalities (38) and (39) are employed to test whether 
Wa. xy) - PISAP, if P<1/2, (40) 
or 
Ha, Xy)- OLSAQ, EOS 1/2, respectively. (41) 
ideally. we would like that (38) always implies (40) and (39) always imphes (41). This wall not always be 
the case because of the approsunations that are used to get (38) and (39). In general however computer 


results indicate that (38) and (39) work quite well in establishing (40) and (41). 


We derived (38) in the following way. Let 2, denote Plax, ) Taking two terms of the Taylor series 
e\pansion of Pebout wy. we lave 


(Pe PLVPY Te ey PY Py. t42) 
: UMaine the expresstons 
| 
: Pra vd [Rta xa] da, vd, Pa xy Rey, {day 
1 
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Wu.) nS. . [lp 202} 
: Maxy? vi, vets »p 2a. 
i tly qed 
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i ge: ry ! 
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1 j x 
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Using this result in (42) 


: oe (P- Py Py = (Py Py WX -xXq) = _ (atl -xy)x-xXy xy, (Xy <a tl). (45) 
Now using an approximation to (37) in (45), namely 
bailing Sone sy tty (46) 
— bole aulayeeyl bas 
we have 
IP- PyliPy = Z — (atl -syle-}-xylley -xy Py, 9 xy Sat. (47) 
Constraining the right hand side to be no greater than A yields (38). 
Inequality (39) is obtained in much the same way. Again 
(0 - Oy WOy = (OVI Qu Xx ~ xy). (48) 
Approximating Qa, by (19),and differentiating (4), we obtain 
lOV/Oyl = (Xy +2 -aiixy tl). xy 2a -2, (49) 
which when used in (48) gives 
IO - Oy 10x = Waxy + 2 alley + I)tly xy. (50) 
Thus approximating (37) as above by (46) and using the result in (SQ) gives 
-Al | : ; 
Oe = = (vy tl eadlyy tT eallvy - vy ly Py, +). Ny aed. (51) 
is . Construming the apht hand side to be no greater than A gives (39). Since the approximation tor Q. given by 
(19), taproves with mereasing x, we used (39) lorxy eta + L rather than ay Oa = 2. 
i : 4. COMPUTER PROGRAM 
; : ‘Phe Haw chart at the end af this wetion on page UP shows in sane detail the baste ansimization and 


How wt the progeam. The boxes marked (A) and (8) reter to particular fonny of the GV +R} whteh ate wed. 
ihus (Ad eters to the algonthn 


Nyey X_L -UMay, PVR EOE, 48 


and @) reters to . 
Xy4y =X, [1 + (Ola, 4,,) - OVRG@ x,,)I. n=0,1,..,° (53) 


P where it is recalled P, G.a are given and R(a, x) is defined by (14). For P< 1/2, (52) is used and when 

P> 1/2 (53) applies. By a single call to the incomplete gamma function ratio routine(PAX 6), [2]. the 
et quantities Pla, x,,), O(a. x,,) and R(a, x,,) are computed to an accuracy of one unit in the twelfth significant 
digit. (There is an option of obtaining 6 or 3 significant digits. instead of 12, if desired.) Ifa < 100, 
R(a, X,,) is always available from (PAX 6). In arder to obtain R(a, x,,) when a 2 100 and 3a < 4x < Sa 
however, an appropriate modification of the original version of PAX 6) was made. 


“¢ The conditions under which (V - R) is terminated were discussed in the previous section. Essentially, 
ua the iterations continue until 


sy iXy ~Xy 4! AY N° (54) 


Then for xy <@ + 1, (38) is checked, or if xy 2 a+ 1 then (39) is checked. Denoting the final output 
value for x by x*, if the applicable inequality, (38) or (39), is satisfied then x* = x,,, if it is not satisfied 
then one more (N - R) iteration is carried out and x* = xy4,- 


Tne provram is presently set with * = 5 X 10-5, By arguments given previously this insures close to 
8 correct significant digits for the output value x*. In the case of the table in Appendix A, A was set to 
id S10, 


“Y A reflection of the very 2oou initial approximations (Ng) for x is the fact that no case was found where 
more than 4 iterations was required, ircluding the 72,000 entries of the table in Appendix A. The methods 
by which x is obtained were discussed in Section 2, 


An estimate of the average computing time per case ior the CDC 6700, with A = 5X 1075, is about 4 
milliseconds. 


A feature of the program permits the user to specify a starting value, As of his own. Thus if the user 
2 desires more significant figures in.x*, he need only make an additiona’ call to the program specifying Xp as 
ge the last valve of. y* obtained, Orit maybe that the user can, in special circumstances, actually supply better 
starting values than those given by the program. 


j As a result of this option tor the user, the possibility ofa poor choice for Xp arises. For ey ample 

i Nii <Q ot xy such that x, <0. or vp such that many iterations of OY ~R) are needed. With this in mind, 

an we a: output ieutian. call i is allotted to inform the user of errors in input or poor choices for 
. tr parncular 


H Ite Yon, a positive integer, indicates the number or iterations by (NV - &). 

‘ Woe then enthera 2 borg <1, Boe $5, (see (and axy flat EP EO, 

: wf 1d. Qo P 20,04 1. In every case an ucceptable value for vy is given, 
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Wede YIP te Ue 88 10-4 ne, Pand G © b- Pare not compatible. Now is found, 


3 
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If - 6 € Y, 20 iterations of (NV -R) have occurred. A poor choice of x9 by the user is indicated. x* 
is given but it may not be very accurate, 


If- 76 ¥.x, from (N ~ RK), using (52) or (53) is negative no value for x" is given. 
It is the user’s responsibility to check location Y. 


The program has been validated for values of a ranging from L0~6 to 106 with Pe [e, 1 -€].(e = 5X 
10-13). The program experienced no difficulty. Also Harter’s table of percentage points [4] were used as 
acheck. A list of the errors found in his table is given in the next section. 


Finally. we remark that it is possible to speed up the program by requiring fewer correct significant 
digits inx*, For example for X= 5 X 1072, the average computing time is reduced from 4 to 2 milliseconds 
with x* pood to 3 or 4 significant digits. A further decrease in computing time for A = 5 X 1072 can be 
achieved by requiring less accuracy in the computation of 2, Q, say 6 significant digits instead of 12. This 
is accomplished by setting a parameter in the PAX 6 subroutine, (see [2] ). Under these circumstances the 
average computing time per case would be about 1.5 milliseconds. 


5. ERRORS IN HARTER'S TABLE 


Harter’s table is a tabulation with $290 entries of 2 in terms of p and P(vly2).. It is contained in {4. 
pp. 137-161]. The input values of v and Pivly?) that he used also make up part of the input for the much 
larger table in Appendix A of this report. 


A careful check of Harter’s results shows that 139 of his values for y? are in error by one unit in the 
6th significant digit which is also the last recorded digit. Certamly such errors detract very little from his 
table. We list the necessary corrections below. The page number in [4] is given first as a herding. The 
first, fourth, and seventh columns sreter to values of p, the second, tifth, and eagth columns refer to values 
af Ay. y2). and the third, sixth and minth columns to the assocuted error, Ifa plus one (+L) occurs in the 
error column (denoted by bh) this means a one should be added in the last place of Harter’s tabulated result 
for y2can case a minus one (- PP oceurs a one should be subteacted. For example, the listing below indicates 
wi error on page 161.2 5 323, Ah?) 8500. Its shown an the error column aya (-1). Therefore Harter’s 
result of yo © 364.848 should be changed to 7 = 304.897. 
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APPENDIX A 
TABLE A 
PERCENTAGE POINTS OF THE CHESQUARE DISTRIBUTION 
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PERCENTAGE POINTS OF Q(ylx*) -- TABLE A 


This appendix contains a table of x2 asa function of vp and Gtvly2). It isa much larger table th nthe 
ones contained in {1, pp. 984-985]. [4, pp. 137-L6L]. There are 72,000 values of x2 tabulated to six 
coriect significant digits with 1200 values of y and 60 values of Q(vh,2) taking the following values: 


B. O.1(0.1)13(0.2)26(0.25)42.5(0.5)5 12 


Q: 0.999999, 0.999995, 0.99999 , 0.999958 0.9999 0.9995, u,.999, 0.9975, 0,995, 0.9925, 0.99, 
0.9875, 0.985, 0.98 (0.008) 0.95, 0.925, 0.9 (0,05) 0.5, 0.475. 0.45 (0.05) 0.10, 0.075, 0.05 
(0.005) 0.01, 0.0125 (0.0025) 0.0025, 0.001, 0.0005, 0.0001, 0.00005, 0.00001, 0.000005, 
0.000001. 


OQ = 99785 and v= 1.9. we have 


The valucs of y2 are expressed in normalized form. For example tor G 
= 233478 Ol = 233478 X 10!. 


y? = 387377 - 02 = 357377 X 1072 or for 0 = 9925, p = 9.9 we have x 


if 
2 


An extremely condensed version of this table is given in Appendix B starting on page 164. 
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APPENDIX B 
TABLE 2B 
PERCENTAGE POINTS OF THE CHI-SQUARE DISTRIBUTION (Condensed) 


Pe . PERCENTAGE POINTS OF Q(v|x2) — CONDENSED TABLE B 


In this appendix, we give 1 condensed version of the table in Appendix A. The valucy of @ and v are 
taken from the table given in [1, pp. 984-985]. There are numerous errors in that table. One of the largest 
occurs for P= 7,0 = 995, They give x2 = 989265 instead of a correct value of x2 = 989256. 


The values of x2 are given to seven correct significant digits in Table B. 
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APPENDIX C 
PROGRAM LISTING ~ GAMINV 


liv. 


GAMINV — PROGRAM LISTING FOR COMPUTING x FROM THE INCOMPLETE 
GAMMA FUNCTION RATIOS AND a 
Purpose 
This subroutine, fora given P, Q and a, is used to find the x value for which 


x a0 
p= i ete har] | e4-ldp a > 0, x20 
0 0 


Calling Sequence 

Call GAMINY (a, x, Xp, 2, Q, TERR) 
where a,x. P, Qare as above xy isan optional input starting value. va must be set to zero if user does not 
want to specity a starting value. IER ts a parameter set by GAMINV. Normally TERR is an integer from 
Ito 19 which denotes the number of iterations the Newton-Raphson method has taken to compute the x 


value. TERR ts stored in ¥ (See page 10). 


IERR > @ when enthera = bora <1} with GMa) 45 eee +ay< 1-4 or Pat, =r 
Pe 0. 2 = 1. dn every case ant acceptuble value ot x is given 


TERR =~ 2 wheng SO. New value ts given, 
ABRK = - d when LP @- FE > $1001, Now value os given. 
fe . VERE =~ 5 when wy SO. Now value is piven, 
TERR = = 6 when the Newtr Rapin method has taken 20 
erations te compute ite V value The x value . 
giten hut ray net he very accurate. 5 a io 
HERR & 7 when usee tras speceticd a peare input Marding value, 


> Asubsucuett ¥, of tha No wt Rapid meets Hand bs 
twgatise, Now value ts ever. 


OMe ae 


Tesimngion peicedute “the tenitine usts tte Newton Rayhson, nicthiod and termites whic m= 20 (which 
Miould faved eseui per whe 


ineey aris Lax 10-8 Seep “aN i9 


” 307: 


we TIARA Mi an tedene 


and one of two other inequalities (38) or (39) are satisified. If (38) or (39) is not satisified*one more 
iteratica is done, 


Accuracy: One can expect about 8 correct significant digits in the output value of x. This can be improved 
by another call to the routine with that x specified as x9. 


Restrictions: a > 0, x9 20 
WP =0 and Q= 1 then x = O is returned. 
If P= Land Q = 0 then x = 10320 js returned [P+ Q - 11) <5°107)3, 


Remark: GAMINV uses 3 subroutines, PAX6, EVA and GAMCQ3. 


Jon 


Program Listing - GAMINV 


SJIWOITINE GAMING € AyXNPL,»S3,PL,NTIT 
OLYENSION C2 Cu) 
DATA ( CACTI pL edpw 2 SelB, oe Dds oC hy olkd / 
DATA RI 271.4142 13562 373107, AL AMS cE @S/ pTALAM/ 1 cEoUs sEPSS Le Em5/,y 
1 CPST Let LaseEPSLSS5 ce EOL 34, I5/9/ 
XNPizl. 
Tfa9Q 
IF € PLeEQeOee AND. QLeEQeL. } VWETURN 
IF ( PE eNEot ee AR eAleNee De ? GO TO 2uit 
XYP121.64 220 
RETURN 
2011 CONTINUE 
TF ¢ 4.67.3. ) GO TO 2:15 
TTs-2 
QZETIRN 
2315 CONTINYE 
IF € ABSCPT+OTaL.) LE-EPS1L ) GO TO 2:47 
Ifrey 
RETURN 
217 SONTINVE 
ITF & S3eGEe*e » GO TO 2d61 
Tte-5 
RETURN 
2061 CONTINUE 
IF € GeNEcde 2 GO TO 251 
XNPL=-ALOG(7I) 
RETURN 
201 CONTINUE 
IF €(€ S3.EQs:. 3 GO TO 3312 
KNS53 
IF ( PI.LE0.5 ) GO TO 3017 
39 TO 3317 
3011 SONTINJE 
T3r1./8 
A>c(PI®A* TS) #9 (73) 
XNZAP/ (1. -AP/(AGi.)) 
T32CXNFAd SCAOL eo) 
XNFLEXKN yas 
IF ( TILTEPS3 } RETURN 
IF ¢( KNeLEede 57 ) GO TO 3915 
BILNzALOG (81) 
ASt4e-lL, 
TLE-BLUNG A3* AL OG(-BiLN) 
X4=-BiL N= (AS) SALOG(T1)-ALOG( (TL42.-AD/ (T1910) » 
XNSAMAX1 ( XN,X4) 
3015 SONTINJE 
IF € PEeGTeOe5 9 GO TO 3317 
3017 CONTINUE 
Irs) 
3024 SONTINUE 
CALL PAX6 ( AgXNoPsQeJ92 ) 
IF ¢€ P.GTeie 2 GO TO 90713 
IF € Pe EQ ee cOReQeEQe6e. ' GO TO 9671 
“KNPLEXNG (1.= (PePII/R) ; 
Tlselfot . 
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Program Listing - GAMINV 


T? ( 17.GEe20 } GO TO 309% 
Teed35(XNFI=XN) 


1F t Ti eG. CALAM®XNPL)) 69 TN 3u34 
Tisfaete 

YF ( GLES(XNPL=i.) 2 30 FO 3326 

aPysbot,. 


TISTLECAPLEKNPLIFABSCA-XNPI OL.) 
TSECTOL AMFXNPL*®XNPS PAPLI A 
IF ( T3eLE.TS 2 REPURN 
G9 TO 3329 
3G25 JONTINUE 
T3TL" (XNPLOL, ADP OXNPLO2.°A) 
TSETALAM® XNP1L *® (XNP 161.) 
iF ( T3eLE.7S 9 RETURN 
3629 SONTINJE 
CALL PAXG (€ Ae XNP1L»PsQedeR ) 
IF € PeGTete 2 GO TH 3971 
IF © 2eEQeds > CALL EVACA,XNPL,R ) 
KXNEXNPY 
KNPLZKNPL* (16> (PePI)/R) 
TTzITes 
IF € I7,.G5€.2: ) SO TO 3,.3 
RETURN 
30271 SONTINUE 
XNEXNPL 
39 fO 3V2t 
T™w2 FONTINUE 
IF (€ 866.456 3 GO TO 5331 
COLL SAMCO3 ( AgHsEPS»TS ) 
3LsQ1°TS 
IF € 31eGE.0.4>5 ) GO TO 3uil 
3t Nz AL0G (31) 
A3zA-k. 
TL e“-31LUN¢43* ALOG(-810N) 
XNZOFLLN@ (AZ ISEALOG(TLISALOG( (TL 42 --ADsS(T1OL1.) ) 
IF € KNeGE.0.67 } GO TO 3316 
T3=1./A 
APzC(PI*A*FTS) **(T3) 
K4SAP7( 1 °APSCA41.)) 
XNZAMAXACXN, XM ) 
3316 SONTINUE 
IF € PleLEese5 } GO TO 3.17 
3317 CONTINUE 
Its) 
3321 CONTINUE 
I= € PeGTeie » GO TO 9U71 
TF € PeEQeOe cOReQeEQeOe } GO TO 9071 
XNPLSXN®O1.¢¢9°QI)/R) 


ITslres. 

IF ( £7-GE.20 » GO TO 3609 
TLSABS(XNPL-XND . 
IF ¢ v1 GT o(ALAMFXNPLD) GO TO 3331 
TisTieTrs Ege 


IF € MoLE(XNPLin-i.) 9 33 FO 3325 
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3325 


3327 


333% 


$331 


5669 


me 


Program Listing ~ GAMINV 


aPLadete 

T3ATL® CAP LEXNPLD FAB SC ASKNPIO@1.)} 
T32CTALAM®XNPL *XNPS PAFLI SA 

IF ¢ TaeLE.TS ) RETURN 

69 fT) $327 

CONTINUE 
T3ISTLOCKNPLOL. “Ad Se CXNPLE2.2A) 
T5sTAlA M® XNPZ*® CXNPLO1 od 

TE ( TILETS ¥ RETURN 
SONTINUS 

CALL PAX6 ( AgXNPL P,Q Jee ) 
TE ¢ PeGTele } GO TO JU71L 

IF € 2%.EQ.06. ) CaLtL EVA ( AgXNPI4R ) 
XNSXNPL 
MNILEKNPL* (Let (Q-Q1) /R) 
ITs(t et 

IF € If, Gete20 2? GO TO 3939 
RETURN 

CONTFINJE 

XNEXNPYL 

GO T3 3321! 

CONTINUE 


P21 


Qz2I 

If ( 31672065 2 THesaatle2-2ALOG(PI) ) 

IF ( JIeLEM.5S ) TsSORT( -2.°ALNG(II) 

KPSTOC( (2621352 34957 L582 TPT 946223 - 21559 B71L'bdFT 
“LleG6L 67202 834965) *Te3esL12 Fs221 Carers 
COCs UOLL7 Oadet 89422-19741 62736 HOSS? AZ2Z239°T 
564069 15977 636349) *Teoe0luS 37656 25662) *T ete) 

If ( PY eL Ee 00s ») XP==-KP 

KOS QzKPeXO 

SLEXP*2AT2 

G2%.265560 66666 660667" (XPSQ-+i. } 

XIC3=2KXPSQSKP 

93ECXP (3 = 7e* XPS (3-FRT2) 

XI eT4=xX FOR*KP 

Ge F~05 eo XPUTH O14.4XPSQ= 32.)/4156 

G52 (9. FXPHTHOEXPH256. MXPC IU SI. *XP)/ (4860 oP RT2) 

Tis2.%A 

RTZ2AzSOART (TL) 

XN805°C TL GL *RI2ZA+G2eSS/2T2AbGe/ TL +G5/(T1 SRT 2A) 

IF € PL.3T6065 } GD TO 5e71 

IF ¢€ 4.676139. 3 GO TO 3017 

IF € KNJGT.C,12 SCAGL009) GO TO 3C17 

IF € heGTewe } GO TO 5469 

AP252A4,5 

L3SINF (APPS) 

TF ¢ PI-GT.C2CL3) 3 GO FO 3C%7 

CONTINUE 

COLL SAMCO3 ( AgH,EPS,TL } 

T3sie/a 

X5a(PIOR* 71) 99 CTS} 

TESA (Lele oOXSI/S (ACL) ) 

MNBCZe9XSI7 (Le OTS ) 
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GQ GO 3917 

SONTINUE 

TF ¢ &,G6€.14.9 9) GO TO $317 

T= ¢ ENGLE (HBA) ) GO TO 3317 

CALL SHMCO3 ( ASH ,EPS,TL D 

31s areri 

BLLNZ 40906 (81) ; 
KVECBLIN| (CL oA PALOGCXND MRL OGIIKXN 62,°A) /SOXNO1.)) 
G) fo 3317 

CINITIWe 

lTs- 

RETURN 

CONTINUE 

Ilse? 

~NOLaXy 

RETURN 

zw) 


Program Listing - GAMINV 


SJRROITINE PAXG & Ay XyANSyDANS,TOPS,& } 

OEWENSION VA6E9),VB020),Z24 030) ,ZB079) 4SA020) $9120) 
9. FOUMAT ( Antia de Aw sE22ehby Sh XE yE22044 ) 
yi FORMAT (45H) TRPUT ERROR ) 

oars 


9U/2.%266 73552 35532 Ee2 43 
PA/2013979 261614 82942 ert 4s 
P276-99n3 83458 o1914 Eog /s 
P37-3.35609 84370 18154 2-2 sy 
4u/72,4505 68758 O3861 EH2 4s 
QL QnL tbe BO5eu 45144 cot 4s 
Q2/1.5062 79753 O4075 Fel 45 
Qivi. 4s 
Cu/3e5045 92610 24162 Er2 / 
QarTa 


174.51 SL 89537 118673 Ce2 ‘> 
027323932 O9107 343 ~HEo?2 i) 
037165296 Sud50 46900 Et2 45 
wefGe3LEZ 22722 20567 =i /s 
Co¢7ecil7? S825G 84559 A+) Cz) 
097503619 S517K 78974 Fol /» 
C7J4.308€ 6857 S$ B2717%O7 Fy 
BL 346/02,9961 7077 V5542 En3s, 
RLPH1L/ 6.9073 U9I06 23251 c-2 / 
OLTh 


Danaus iw we DBT se Ge MK 


1 ALPH27 =f. 2695 05435 39097 Feil Zy 

2 AL7H3 79207066 13046 C9648 Zl sy 

3 GLPMG/926 2319 24597 S4id5 E-2 /y 

& 35001.5620 923.5 26468 Ev? a) 

> 82971/1.9150 63261 & 7830 Eok a) 

4 829F2/13.0516 75107 00793 ‘ys 

7 8293/4,9973 32010 17135 49 

6 SETa/i. a) 

9 02735.30495 94659 50952 Ee? / 

Dora 

& 347729095 09253 27898 Fe2 /s 

2 3279035135 40946 50610 fre 4s 

3 337363699 C2644 65631 Ee2 4s 

4 007207758 Slone? 03988 202 4 

& 0577-7000 15293 52295 cit 4s 

6 95/1.2782 7273% Go29u Fok ‘4s 

7 O7/ i. 4 
c SST28SP(OP3415./0P139398 Bf. ) 
c oSU38SP (OP Joi73S./OP 29155 3920, ) 
c 25 e3SP (07411360 8261-/0P5L791 SB7d— O. ) 
c oST55SP (O0P75652 u2S330/0P 76364 1680f 9 LE. ) 

vata 

1 C5154. 96539 65736 46935F-3 fy 

2 5504/4, 22422 66000 SECI2E-2 4, 

3 C3084 01 79H7 22737 TISQSES2 7, 

& GETZ e243G4 1469 7HLbbEHe 4, 

5S Gh4955/1,7724 53850 90552 4s 

& QE2PT/2.506% 26274 63199 / 

FyMTFt/4 772s COSI 90552 


/ 
DOS. alOL 5 rE 157 p mE S/3 Le 7T sg ALPHAS/ 0.8927, 


{73 


LACE 


1107 


1155 


1161 


137% 


FwN ee 


ALPHAI/ 7. 3488/,ACO/ in ‘5 
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21731-0007) X05/17.007 9CEP/E 43025 35092 994605/ 
9UK141.0366/ ,UX3/2.9205/ 
9°T24 14142 13562 S7310EO1/,A3A9/1006/ 


pALPHAS /3, 290534 9X50 90 7p IX5/ 20 CBT yp ACCS/5 EM Gs 


LSAavsa 
IND3) 


ANS BQ, 


IF € LOP3.NE.. ) 
A-PHA= ALP HAL 
XUEKGL 

ACC=ACC1 

Jus IMS 

GO TO 411411 
CONTINUE 

IF ( ITOP3.NE.1L ) 
hi F4hz ALPHA 
AsC=4C53 

XG2xo3 

UX sUX3 

39 FO Lids 
CONTINJE 

ALPHA ALPHAS 

ACS ZACCS 

SE? On) 

UKEIXS 

CONTANUE 


GO TO 1105 


GO TO 1147 


IF € Meh Eodod cOReXetLToOeu ? 


IF € KeGTeNeu 
CINTINJE 

30 Fo 1171 
ANSz2. 

RETJIRN 
CINTINUE 

IF € AeLT.15~ ? 
REAsSIIT CA) 


GO TO 1151 


GO TO 1331 


IF € ACOoLT-ACC? ) GOTO 1155 
TF ¢ AGE.X ) GO TO 11518 
FIstertes (94) CALPHAS (3e°RTA) 


FPS-F PEF POF P® ASX 
IF € 7P.GE00.9 2 
CONTINJE 

TF € Sel TeASA2 ) 


GO TO 1191 
GO fo 1331 


GO T9 1131 


LE © Ce 7S PAV GT AXeORs CL 25°AVOLTAX «GO TO 1331 


GO TO S611 
sINTINUE 
CINTINVE 


FMslerLes (90% A) MALPHAS(3e°RTA) 


FM sao F MO FHOF MFA OK 
IF @ FHL Eee 

GO TO 1158 

CONT INJSE 

ANS#3. 

GQA%NS71. 

RETITN 


GO FO i821 
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219%) SONTINUE 
131100 CONTINUE 
ANS21. 
QANSz3. 
RETURN 
1331 CONTINJE 
I= ¢ ALEK) GO TO ist 
1341 CONTINUE 
‘ 1F ( AeGOe30e ) GO TO 1323 
ALSX=8L96(x) ae 
TizeoXed SALGX 
CALL SAMCO3 (€ AgtiyACC,GANS ) 
Q2EKP(T1L) /GANS 
90 TO 1349 
(2363 CONTINUE 
AIS VEL ./ (AFA) 
TLeA~KOL*ALOG( M/A) C6 CALSQ01.63333 39333 33333) *AISQ44.656 65660 
1 666667) *AISQ-140.)/( 1080.°%A) 
RERTAEXP (7199239994 22804 91433 
1369 CINTINUE 
T= € beGT eX) GO TO 135% 
1356 SONTINUE 
IF € AeLTeod- ) GO TO 1453 
IF ¢ A&.GE.X ) GO TO 3341 
90 70 357 § 
1351 GINTINGE 
IF © WLE ole 2*AS( Lem X/(A41.))*ACO)) GO TO 1171 
GI TO 135% 
deGL SONTINUE 
ITF (€ X.GEeX%d 2 GO TO L361 
TAdJAz2.°A 
JEIntT (TWO A) 
TIzASSCFL CAT ( J) -THOA) 
IF € P21 .6Te00)60 TO 1341 
L23/2 
AFshe-*LOATI(T) 
TF ( AF GTeUe » GO TO 601t 
a2 79 6071 
4459 = SONFINUE 
TF © Koh T oh 25) GO TO 1461 
TEs (2 0¢X-AN/ 61 6 OX) 
TE ( Re lLEe (ACOPX*T4*%.293G60 TO 1311 
IF € ReGT ol L029X*T1) 2 «30 TO 3311 
69 70 7911 
1461 CONTINUE 
INJ©1 
20 9 3042 


CG 
C THIS TS PART F 
Cc 


2013 CONTINUE 
ACG7s,5°ACC 
F7sA°AlGX 
CEEs2, : ee! 
TST7/SEE is 
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SUszT 

293% CONTINJE 
COlzSrFed. 
TECTOT7/CEE 
SUMSS IM of 
IF ¢( Q@3S0T).GT.ACS7 ) GO TO 2031 
T2377" ( 3. #SuM) 
TKeex 
SEEL. 
AISKBTK/ (ACCEL) 

; SIMZA JK 

2u71 SONTINJE 
TELECZ EOL. 
TKE-(X*7TK)/CEE 
a JK2TK/ (ACCEE) 
SUMSSUMCA JK 
IF € &B8S (ASK) eGEe-ACC* SUM GO TO 2571 
T3e-ASSUM 
Tish 
TetILeor2 
T5sfiet2 
AANSSPS-TGel3* (TH+TS ) -T5 
ANS21.-QANS 
RETIZA 

3UC& SONTINUE 
IF € XeLECEP 9 GO TO 3011 
IF ( Ket TeoXO } GU TO 7012 
IF ( A.LT.2. ? GO TO 3006 
T™ 6 (ROCXOL od doh Eo FACOPXK® (2.eX2A)*®.2))9GO TO 1311 
GO TO 90131 

30CG CONTINJE 
IF © Rel Eo Ce. 2° ACO*(XOL.24))) GO TO 13711 
GO fT) 9011 


THIS IS PART A 


anan 


3012 CONTINUE 
HAVACCSACC®.5 
ROVASKR/A 
SHPLste 
SUMZ5NP1 
ClEsa 

3031 CONTINJE 
CEESCIES1, 
Tisx/CEE 
SHPLSSNPA CTL 


Sttmel ita ee Me o 
wer wwe wit @ 


IF € TheGTe 665666 66666 66667 ) GO TO 3031 
3041 CONTINUE 

IF ¢ SNPLLESHAVACC 2 30 TO 305t 

CEESCEES1, 

SNPS (SNP °K) / CEE 

SUASSUMOSNPL 

GO TO 3041 
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5371 


$362 


5369 
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GONTLNUE 

ANS ZROV A® SUM 

QANS#1.-ANS 

IF € [NOeEQe1eAND.ANS.G7.0.9 } GO TO 201% 
RETURN 


THIS Is PART E12 


CONTINUE 

ATSAsi./ (ASA) 

Tish XeA® ALOG( XA) © (0 (AISQ21.23333 33333 33335) *AISQ+4.666 66666 
666657) * ATSQ-140.)/(1680.° A) 

RERTASEXP(T1)%.3909% 22804 01433 

ASAVSA 

ArAel. 

AL PHY= 3.*SQRT (Ae2./3.) 

ALPINV=1. /AL PHY 

A2THRG= At .66666 60666 66667 

DELTA=(X—A)-.666606 66666 66667 

SZSCELTA/A2THRD tt 

ZEEXP (633333 33338 FSITI*ALOG(X/ACTHRO)) 

SOZCAL PHY SZ) 7 ((Ze1.) *Z¢1.) 

XSAv=EK 

X2(S0°SC) 7/2. 

QIKZABS (SUI/RT2? 

TSZEXP(-X) 

IF ¢( Xe67.0.25 3} GY TO 5361 

TLE 0CP3*X oP 2)*XePL) OX OPO 

T38((Q3%X0Q2)°XeQ1)*X+Q0 

OSL sRIX*CTL/T3) 

SERF81.°D0EL 

30 TO 5364 

CIITINUE 

TF ¢ XeGTei6.) GO TO 5371 

TASCIEATXF(CLERTX# (CZEREKY (CIEVTX® (CHEAT X4 (CH ORT X# (CH C7ERTK) ID) 


’ 
TIZDOCRIX*(DLERIX® COZ 4RTX® COS+ERTX4 (CDGORTK* (OS HRT X® (DGeRTX9G7 IID DD 
’ . 


SERFsTist3 

JEL21.+CE RF 

GO fd 5361 

CONTINJE 

TsLosK 

TA SAL PADS T® (ALPHIET® (CALPH26 T* (ALPHS¢ TAL PHS)? 
TISBET (eT (BETLIOT® CBE T2eT*(BETS¢T*BETG))) 
CERFa(Ze/RIXIS(L/RT PLlOTis CT3*X)) 
OELSL.-CERF 

CONTINJE 

ZAjsoL. 

IF ( MeLE 0625 » ZBOB(.5*SERFFRT2PID/ TS 
TF € KeGT 09025 ) ZBO8.5*CERF*RT2PT 

IF ¢@€ SO,GE.G- » GO FO 53869 

ZAQBWLe 

CONTINUE 

XsxSAy 
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230102735 
ZAC4) 2ZA0 
$352=5)°S) 


rrecuszaasisa COPY AVAL ARE TR Png OES BOT 


02 5391 L=255 PEFR iH Fi outs meray 
Ke2°(L-1) Chinett i Wit nce Seg 
T72T7*S5S0 

Z3CLV FHT? 


CF LOATIKI=1,9°23 (Laid 
ZACLIZ<T7 *SO*FLOAT (KD ® ZA (L=1) 
SONTINUC 
¥3¢1) 225¢1) 
V302)2 ©,63333 33333 333336-1°25(3) 
V401) =.66666 66660 66n07E-1L ®ZA(3) 
$301) =v 9¢1) 
APHSISAL FINVE ALPINV 
ALP AC 3=4L PHSQ* ALPINY 
$342) =¥B( 2) *ALPHSQ 
SACL) 2V401) *ALPHC3 
1155S83¢4) 6832) 
Ti725401) 
V5C39S = 0629060 -2.9625° 23059) 4.55555 55555 55556EW1 
VAC2) 2647619 IG7EL JL47BEW1 O2ZA6 4) 655555 55555 S5550E-2°ZA(5) 
T2328.7HSQ 
1523 
MF s) 
CONTINUE 
T7e15-1 
T23zT23°ALPHSO 
$3015) 3v8 (15) *T23 
$4017) sVACT7)*7T234AL PING 
T1i52T15¢S 3615) 
TAT STZ eS ACT7) 
IF € ASSCSACIS) DV LTeTOLecANDABS(SACT7))-LT-TOL 2 GO TO 5431 
1322°15 
I5sI5e1 
T7e8T7°SCS9 
T23586.°FL CAT (15) 7, 
ZECP 3 s-T7OT2]5°7Z94T $-1) 
24013) oT 7950 (T2501 69%°2601 3-1) 
t7s17°S¢S9 
Z3CT301)2-77407T25¢2.) *23(23) 
ZOCTS 12 z-T7*SOei hese 3.) eZAC1L3) 
Mi aM™ of 
GD TO © F397, 5399 954019 Se 55 56059 5407p SCI SH Lip SHlL3e5HL5ySHSL) aM 


1 
CONTINJE 


VdCOd t= (2805) = 6 1046s eel Hhebe®Z8(6P+eZ51N8 LHBLe BL4EALES2°2307 
dd) *.41666 60666 66457E°1 
VACSDSEZA MS) 626714 28571 4265747A(6)4.00025972A(79)9%. 37037 03703 
7) 370E-1 
G2 TO 5395 
CONTINUE 
V305) 5-(2566) 526379 145862 OLC58*ZA(7)be11342 59259 25926E=1* 2919 


178 


Sk 


™N = 


Ne 


S43 


ft 


nme 


5405 


wre 


we wn 


5407 


wre 


we 


bal 


oeunpr 


oe wh = 


Sell 
2 


a ae 
ef ! o@ i 4! t A 
Pipes te rary ‘ J r “fs 
i : 881 gia Le : 
Chik if i os dies 
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Y-.60281 63500 2669LE94%2809))9. 355332 S3333 SSSSSF HL 

VOCE 65GIOS BDIOIY SOSISES1*ZA(0)-.A5097 10176 3668s E292 Al7) 

652536 bo .67 OL 9GUT@S9Z Alb Heh SOC G1152 26337 cC-54%75 (9) 

$9 19 5395 

CONTINUE 

¥306) 2= (29072 631305 85987 5263497818) %.236C6 21693 12169E91%2319 

92067 .67 390123 HSE7IE@3*ZH(1 26617050 52716 0493 8F45°23011))* 
027777 77777 F777 aEMW1 
VA05) 2.25641 32564 10256691 °ZA(7)-,357892 146960 S2601LE-2°ZA(35)* 
070924 L0225 Te IE OE-3* ZAI. 2156F 19917 93063E-4%2Z4 (19) 
Z740119%.15395 91966 72154E-6 

93 TO 5395 

SONTINUE 

VACZ)F¢.2 2222 22222 22222c-1%ZA(8s 2.68026 17135 Wuv7Emec*ZAC9)e 
eL0129 13268 CBC ue E-2*ZACLi I MeG 2578 29446 FICHBE4*ZA(11)% 
2624637 74167 L3o2GE-6*2Z4612) "622326 53177 B6I2ZSEmo*ZA(13) 

¥B(7)2-,.23809 52360 9523862°1°Z2306)+.88185 92672 75934E=2°2B19)- 
09Cbe2 29253 23279E-397H8(10) 6.31511 18214 77562E-4*ZB(11)- 
032596 73639 HdSC7TESG*Z3(12) teb0513 60737 22756E-6*29(135) 

69 T9 5395 

CONTINUE 

VEA07) 6.19007 Be313 725KU9ER1*ZA(9) —. O6G3E CI35* 26792E-2°ZA LI) + 
eALFG0L G7u27 C3IS7eEm2?>ZA(11I 6566055 77635 7IGIOEM4*ZA(12)+ 
015617 53579 Bde23E259 724113) 2.12799 92297 SSILTE M7 *ZA(1L4)6 
731309 C719 48566E-19%Z28 (15) 

V¥308)=-.25833 33333 33533E-1°23(9) 66875464 99954 2° 356E-2°*Z3(1.)- 
011592 44992 7I2S7EZ6LZIC11D + e5433d GLES2 2175 9F-4*Z3(12)- 
015360 16990 95d43E-5%28(13)+,63258 50678 B294HEn3* CB 14) = 
055373 364275 SOs00E 11°28 (19) 

G) TO 5395 

CONTINUE err 

VACA E61 7563 65904 912262- 1928110) ~ 064985 16229 139KLEW2°ZA (11d + 
o13493 69333 ds0o72ZE@2°ZACL Ee) .92uH1 BO4C7 TOL27EMG*ZA(1 3)+ 
029398 23607 SSSTHEMS®ZA( 14) -0427"7 63569 SIN7SEW7¥ZA(1L5)6 
24541 46564 68959E- 992A (1b)-.56915 56190 33935E-12%24 (17) 

V309)27.26518 SLU5L B5155E91L9Z2(10)¢ 065966 42554 I7564E"2°23t11Lb- 
012740 20646 36950E 29723612) 667923% 62192 EBISZE-4*ZB(1 5) - 
6219863 32483 9>20~E 5925 (14) 6025640 Gohl L2548E-7*ZA(15)- 
016526 266535 56316 £9923 (16) 6.5760" 56531 7TAN2CEL4*2Z3 (17) 

62 10 5395 : 

CONTINUCL 

VACOD= 4.15373 Ulad7 Suld9—o19ZA(11) 9662621 SOL5S 67318 -2°7A CL 2) 

olbO52 16935 OLeTFE—2*ZACAL SIN eL1757 92069 SAILLEHS* Zh (14)t 
047065 7699) 799S4EM59ZAC15) 2.95202 29662 GBIZHE-7*2ZA(16)t 
092392 73029 CHdSLEW3ZACL7I Me 3H4H52 25970 SOSGGELAI*ZA (LODE 
030053 71021 16680E 9149248119) 

¥3016) 5-21666h coobo 606672 -1°23(11) 6.63928 J1741 7322260297292) 

016682 22299 SeuS3Ene292Z3(13) 461485 97493 EH735E-S*7ZB0 14> 
037756 944538 863969E 75920115) +.55759 96057 S7A53E-7*ZB(10)- 
o2b453 EQ739 ObL2VE<9Z3CL7) HL SHL2 S512 LE9IIGEMLI ZS (Ldd- 
093497 93064 941676915923 (19) 

G) 10 5395 

CINTI WE 

VAC10) 86614492 75S02 SLOBHE94%ZA012)9.8053m 46095 24176F -2e2Alls) 

$e15625 JO7OL E2S27PEM29ZAC Law eL 4249 B27H0 VIIDIESS*ZACI5)e 


¢ 


179 


aan 


Se 31 


Gules 


5u5d 


i 
2 
3 
6 
) 


037347 bn 263 
eo2euu7 91370 
ewb225 32193 
v3adtbd=-,15252% 
015164 52345 
or7 Le? Wel? 
01538, 56944 
015549 66311 
60 TT) 5395 
CONTINUE 
VACLLP 6.13333 
@oloeite2 W535 
031269 C9027 
032343 42766 
26317 C716% 
0239671 - 7203 
¥3¢12) e-.13888 
=210039 76493 
079380 ai0ii 
035925 GEbL 
212351 §23915 
033717 12666 
Gd Td 5395 
SONTINLE 
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92697 c-5*Zh(Lodeo1L 7514 13319 SI7TLSE 6% ZA(1794 
HLOLZE99ZAC1L6) 010288 16509 G4713E“109ZA (1906 
27500E$13°ZA (2ude. 35605 28723 32111 E-16%2ZA (21) 
5L525 151525 21°230.2) 2.62685 477395 418286 -2%23(13) 
S28S5E$2S9ZB Clade. 1L3IOLG 45605 38248E -3°728115)- 
O25e9E=5°Z8( Lode. 13158 76549 27232F -69ZB8(17)- 

7437 3E-d923(16)%.8213F O7414 60955E-119%728 (19) - 

67 6O5E$139729 6209 %.466506 60904 15139€-179723 (21) 


33333 SISSTE ML ZAC1 332078226 73159 Y79IKWSE“25ZA(14) 
SeLBIEm29ZAC 1S) el 5641 22539 S3S5S30E“3% ZA(1E)6 
5U3Y7E$392ZA617) 2.28427 17670 98726E-6*ZA(1 Be 
29405E-89Z0 (19) 649118 23691 34031E 10724 (20) + 
7TO~BLE“129ZR021)2251251 26773 74051 E-15*%ZA(22)¢ 
454266 -139°Z24 (025) 

883865 6346892219 29(13) 2.79379 48362 75189E-29759(14) 
89170 b-2%728 615) 0154659 76695 11729E<3*ZB(16)- 
280b0E-5*28(17)6.22641 75615 827386 -6%29(18)- 
66,41LEF9°3°23619) 6.29351 25130 22109-10728 (20) - 
67799€-12%Z8(621)+.15785 01067 339035-15*%2B(22)- 
17530E-13928 (23) 


VO(L22=¢,12345 B79UL 236575 -42*ZA(14) -275962 C3771 T24B0EH249ZACL 
5d %617048 20355 9010 3E-2% 2ZA(16) -018910 05799 70BSSE-S3*ZACL7)& 


011572 9332. 
031776 186615 
094)15 10565 
050384 30725 
03023) 2.12829 
= eib74&7 26884 
014583 54179 
056844 79270 
067532 7376. 
el42u1 618646 
G) TI 35395 
Ss ONTINJE 
Tsi./8 


SoH hGEa*ZA (18) 042273 HH701 A79KBE-6*ZA(1 IF 
L4S62E=5"ZA0 2b) 611653 40493 25387E-9*ZA(A1)¢ 
BIsob1cZh (22)-454318 52910 B5580E-14%ZA (23)% 
L400 dG EM179Z0(Z4)=-.22478 06537 4SO20E-20*ZA (25) 
DL282 O5120E 19239614) +677008 28747 69910F-29728 (15) 
27663 E2423 (LEI O.L77IL 99594 B5527E3*29 (179° 

0550 2E-492B(16) 6234937 60792 262625 -64%720(19)- 
TLSB4EMd* ZB C2U) Fe774H94 GSO3CO 60715E~-10%25 (21) = 
267156912928 (22) 4.13508 28991 S8372E-14%7Z3 (23) - 
B9S2SF “17978 624) + e23444 67226 51062F-21*%28 (25) 


CSTaC( CCC CSIS*T CST 4) *TOSST 3) FT eCST ld FT .11959 87054 S32: 99E°1) FT 
o27777 77777 F777B Ent) *Tol. 
POPLXSICST*CTI56017)*T5)7K5 2P1 


AzASAd 

IF © Se ok Tede 
QAaNSsPAPiX 
ANSsi e-JANS 
GO TO 3451 
CONTINJE 
ANSsPAPLX 
JONSsL.-ANS 
CONTINU: 
RETURN 


’ 


GO TO Sei 


415 IS PART 61 


6011 


6034 


6071 


6091 


6101 


‘141 


Oana 


7041 


7U38 
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CONTINUE 

NzQg 

RTKsSQRT(X) 
TSSEXP(-X) 

RNET5S/ (REX*GANPTS) 


THIS IS PART SL SUPPLEMENT 


IF ( €eGTet16.) GO TO 0331 
TLECUORTX*(CLERTX® (CZeRTKX* (CIRIN® (CUERTX* (CHORTX*® (CO*C7TFRIXI DID) 


) ‘ 

TIZOOERTX* (OL CRIX* COZHRIXK® COS*#RIX® (OGERTX*® (OSERTX® (D64RTX*0P) 109) 
) 

OELECTS*T1) /T3 

G23 FO 6291 

CONTINUE 

T=h0/X 

TIL SALPHOOCT# (CAL PHLeET® (ALPHZ+ T*® (ALPH3*T*AL PHS) )) 

TISZJETOCT*(ZETLeT* (CSET2¢F* (BETS+T* SETS)?» 

OEL=s¢ 15 SUTXIFACLe/RIPI *¢Tis (T39X)) 

60 739 5091 

CONTINJE 

Net 

RNSEXP(-X) 

DZusRN 

CONTINUE 

PNEXEDEL 

CESZAFOFLCAT (ND o1. 

CONTINUE 

IF ¢€ NeEQeI ) GO TO E4114 

NeVeot 

CEEZSE Cele 

QNS(X*RN) /CEE 

PNOXZPNCX CRN 

60 TO 6193 

sONTINVUE 

R=K*RN 

QANS=£9NCX 

ANS#1.,.<-QANS 

RETURN 


TAIS IS PART C 


CONTINJE 

A2NML21. 

A2NS1. 

B2NML SX 

B2NsXe1.eA 

aNsit. 

CONTINUE 
A2NNLEX*® A2N¢C ANS AZNML 
32NN1LSX°B 2NOANTBZNML 
ANDO SA2NMi /B2NN1 
ANSANO1. 

SA2NZAN-A 
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A2NZA2NM1L CSA2NFAQN 

B2N2B2N41 +¢SA2N*B2N 

AN} 2A2N/52N 

IF (€ ABS CAND AMG) GE. CACC*ANG) 9 GO TO 7032 
QGNS=2FANO 

ANS=1.-QANS 

RETURN 


THIS IS PART O 


ono 


QuiL SONTINUE 
SNP i=, 
SUMZSNPY 
CEE=A 

9031 SONTINUE 
CEESCicel. 
SNPL3(SNPA®CEED /X 
SUNSSJM+SNPL 
I= ( A3SCSNPL).GE.ACC ) GO TO 9031 
QaNS=(R® SUM) /X 
ANS21. -QANS 
RETURN 
END 
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3011 


3931 


3651 


3071 
3093 
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SU3SROUTINE GAMCO3 ( AygHy, ac ,Tt ) 

OLSENSION C(30) 

OATA ( CLI) ,Is3,18 ) / 
#.57721 56649 013533 cod 
055587 60715 232546 Eeg 
meh 2002 63583 46952 cmt 
#616653 86113 82291 E+0 
2e62197 73455 55u43 Enk 
295219 71527 67697 E-2 
o72169 ©3246 66311 E-2 
“211651 67591 05907 E-2, 
-e02152% 16741 14951E-3, 
#012605 02623 68116 €-3, 
220134 854786 47852 T-% 
#012504 93492 16267 Fes 
#011330 27231 98176 E-5 
220563 384616 97761 E-s 
#535160 9510% «8442 Eo8 
#250020 0764% 46922 F-8 

A25e.5%AC 

T2 INF (a) 

A= she-FLOAT(I) 

IF ¢ AFCO - 2 «GO TO 34391 

CESzAFHe- 4. 

VF C GF T2065 2 CEERAF 

PHISClE ; 

<z3 

Ands0 (4) CEE 

SUMZANK 

CONTINUE 

CEEZSZE*PHI 

KeKel 

ANK8C (K) ®CEc 

SUSSSJMe ANK 

T= ( BIS CANK) cL Te GACCEARS(CSUM) > 9) GC TO 3u31 

IF ( &.EQ.15 ) GO TO 3031 

GO TO 3011 

CONTINUE 

HsSJ" 

IF ¢€ 47.GE.0.5 ) Ha(1.°AeSUM) ZA 

Tish. SUM 

IF € AFL T0005 ) TLISAF*T1 

Tiste/T1 F 

IF ¢€ 1.€0.0 ) RETURN 

AJSAF oi. 

Js) 

CONTINUE 

sujet 

AJSAJOL. 

Casas*Ti 

IF ( Joh Tol » GO TO 3uUSL He 

CONTINUE 

GO 7) 3151 

sOVS INUVE 

PEP | 


ww @ewewe es 


sew ewe 
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Tistie 

IF ¢( IT.€Qei } RETURN 

Isle 
3131 CONTINJE 

JeJeot 

TLsFLOATCJD9TL 

IF (€ Job Tol } GO TO 3138 
3151 CONTINUE 

RETURN 

EN) 
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SJFROUTINE EVA ( AyXyR } 

DATA ACC/5.E-137 

IF € QceGEe30e 2 GO TO 1303 

Ti sex +A%A 0G (x) 

CALL GAMCO3 ( AsHsACC,GANS ? 

RIEKPC TL) /GANS 

RETURN 

CONTINUE 

QTazSQRT (AD 

ALSAzLe/ CAPA) 

Ti zhe XAFS ALOG(X/AD © CC CALSQ"105333 JFTIZ Z%33S33)AISNCHOH6 66466 
666657) * AISQ=1 40.) £61680.°A) 

RERTA*EXP(T1)*.3989% 22806 01433 

RETURN 

EN) 
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